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ABSTRACT 

We report on a multiwavelength campaign on the TeV 7-ray blazar Markarian 
(Mrk) 421 performed during December 2002 and January 2003. These target 
of opportunity observations were initiated by the detection of X-ray and TeV 
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7-ray flares with the All Sky Monitor (ASM) on board the Rossi X-ray Timing 
Explorer (RXTE) and the 10 m Whipple 7-ray telescope. The campaign included 
observational coverage in the radio (University of Michigan Radio Astronomy 
Observatory), optical (Boltwood, La Palma KVA 0.6m, WIYN 0.9m), X-ray 
(RXTE pointed telescopes), and TeV 7-ray (Whipple and HEGRA) bands. 

At TeV energies, the observations revealed several flares at intermediate flux 
levels, peaking between 1 and 1.5 times the flux from the Crab Nebula. While 
the time averaged spectrum can be fitted with a single power law of photon in- 
dex T = 2.8 from dN 1 /dE oc E~ r , we find some evidence for spectral variability. 
Confirming earlier results, the campaign reveals a rather loose correlation be- 
tween the X-ray and TeV 7-ray fluxes. In one case, a very strong X-ray flare is 
not accompanied by a comparable TeV 7-ray flare. Although the source flux was 
variable in the optical and radio bands, the sparse sampling of the optical and 
radio light curves does not allow us to study the correlation properties in detail. 

We present a simple analysis of the data with a synchrotron-self Compton 
model, emphasizing that models with very high Doppler factors and low magnetic 
fields can describe the data. 

Subject headings: BL Lacertae objects: individual (Mrk421), galaxies: jets, 
gamma rays: observations, radiation mechanisms: non-thermal, X-rays: indi- 
vidual (Mrk421) 

1. Introduction 

The space-borne EGRET (Energetic Gamma Ray Experiment Telescope) detector on 
board the Compton Gamma-Ray Observatory discovered strong MeV and GeV 7-ray emission 
from 66 blazars, mainly from Flat Spectrum Radio Quasars and Unidentified Flat Spectrum 
Radio Sources (Hartman et al. 1992, 1999). Ground-based Cherenkov telescopes discovered 
TeV 7-ray emission from seven blazars, five of which were not detected by EGRET (Horan 
& Weekes 2004; Aharonian et al. 2005). Although 7-ray emission from blazars have been 
studied for more than a decade now, it is still unclear where and how the emission originates. 
According to the most common paradigm, the emission originates close to a mass-accreting 
supermassive black hole, in a relativistically moving collimated plasma outflow (jet) that 
is aligned with the line of sight to within a few degrees. The relativistic Doppler effect 
can explain the intensity of the blazar emission, and its rapid variability at X-rays and 7- 
ray energies on hour time scales: for emission originating from synchrotron or synchrotron 
self-Compton (SSC) models, the apparent luminosity increases approximately as the fourth 
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power of the relativistic Doppler factor 1 <5j, and the observed flux variability timescale is 
inversely proportional to <5j. 

Blazars are powerful sources across the electromagnetic spectrum. Typical spectral en- 
ergy distributions (SEDs) for the high energy peaked TeV blazars show two broad peaks, one 
at infrared to X-ray energies and the other at X-ray to 7-ray energies. The low-energy peak 
is commonly believed to originate as synchrotron emission from a population of relativistic 
electrons gyrating in the magnetic field of the jet plasma. The origin of the high-energy peak 
is unproven. The commonly adopted and best studied models assume that the 7-rays are 
produced in inverse Compton processes by the same electrons that emit the synchrotron ra- 
diation at longer wavelengths (for a recent review of observations and models, see Krawczyn- 
ski (2004)). In so-called hadronic models, 7-rays are emitted as synchrotron radiation of 
extremely energetic protons (Aharonian 2000; Miicke et al. 2003), as inverse Compton and 
synchrotron emission from a Proton Induced Cascade (PIC) (Mannheim 1998), or from 
7t° — > 77 decays following the interaction of high energy protons with some target material 
(Pohl & Schlickeiser 2000). Recent reviews on observations of blazars with TeV emission 
and models developed to describe the data can be found in various review articles and books 
(Krawczynski 2004, 2005; Tavecchio 2004; Weekes 2003; Aharonian 2004). Reviews fo- 
cussing on observations and models of sources with MeV/GeV emission can be found in 
Sikora & Madejski (2001); Coppi (1999). Broader overviews of the field of TeV 7-ray 
astronomy are given in Buckley et al. (2001); Ong (2003); Weekes (2003); Aharonian 
(2004). 

Multiwavelength observations are key for understanding the blazar phenomenon. The 
acquisition of good multiwavelength data sets has encountered substantial difficulties as the 
sensitivities of current TeV observatories require flares for sampling the TeV light curves on a 
time scale of hours. Some sources were observed with excellent multiwavelength coverage but 
during relatively unspectacular quiescent phases; in other cases, the sources were flaring, but 
the fluxes were only poorly sampled in frequency space and in time. The most remarkable 
result from the multiwavelength campaigns is that there is good evidence for a correlation 
between the X-ray fluxes and the TeV 7-ray fluxes for the two sources Mrk 421 (Buckley et 
al. 1996; Takahashi et al. 1996, 2000; Blazejowski et al. 2005) and Mrk 501 (Djannati-Atai 
et al. 1999; Sambruna et al. 2000; Krawczynski et al. 2002). 

In this paper we present results from a multiwavelength campaign on the TeV blazar 



1 The relativistic Doppler factor is given by <5j = [r(l — f3 cos^)] -1 with T, the bulk Lorentz factor of the 
jet plasma, 0, the angle between jet axis and the line of sight, and /3, the plasma velocity, in units of the 
speed of light. 



- 5 - 



Markarian (Mrk) 421. The source is a nearby (z = 0.031) high energy-peaked BL Lac object, 
and was the first extragalactic source detected in the TeV 7-ray band (Punch et al. 1992). 
In November 2002, 7-ray observations with the Whipple 10 m telescope revealed several 
Mrk 421 flares with fluxes exceeding three times the steady flux from the Crab Nebula. The 
All Sky Monitor instrument aboard RXTE also showed extremely strong 2-12 keV fluxes 
reaching 100 milli-Crab. Collectively these triggered a coordinated campaign. We invoked 
radio, optical, and X-Ray (RXTE) observations to commence as soon as the waning Moon 
would allow the Cherenkov telescope to take data once more. Although the X-ray and TeV 
7-ray fluxes had decreased substantially when the campaign started on December 4th, we 
acquired a data set with a high signal-to-noise-ratio X-ray light curve and X-ray energy 
spectra, and good signal-to-noise-ratio TeV light curves and TeV 7-ray energy spectra. The 
combined data allowed us to study the X-ray/TeV 7-ray flux correlation over several weeks. 
Following our previous study (Blazejowski et al. 2005), this is the second campaign that 
measures the X-ray/TeV 7-ray flux correlation over several weekes. For the first time, we 
use here simulated lightcurves to address the statistical significance of the X-ray/TeV 7-ray 
flux correlation and to constrain the time lag between the two light curves. Simulations are 
necessary as subsequent data points in the light curve are not independent of each other 
(see e.g., the discussion by Edelson et al. (1988)). The rest of the paper is organized as 
follows. After describing the data sets in Section 2, we explain the method that we used to 
reconstruct TeV 7-ray energy spectra from the Whipple data, and give the results obtained 
with the method in Section 3. Subsequently, we present the results of the campaign in 
Section 4, and conclude with a summary and a discussion in Section 5. If not mentioned 
otherwise, errors are quoted on the level of one standard deviation, and upper limits are 
given on 90% confidence level. 



2. Observations and Standard Data Reduction 

Radio Observations 

We used the University of Michigan equatorially mounted 26-meter paraboloid in its au- 
tomatic observation mode to observe Mrk 421 at 4.8 GHz, 8.0 GHz and 14.5 GHz between 
December 3, 2002 and January 10, 2003. Both linear polarization and total flux density mea- 
surements were made, but only the total flux density measurements are reported here. Each 
observation consisted of a series of five-minute ON-OFF type measurements over a 40-minute 
time period interleaved with observations of other program sources. Observations were re- 
stricted to within three hours of the local meridian to minimize instrumental errors. The 
telescope pointing corrections are interpolated from position scans through sources stronger 
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than 2-Janskys. The flux density measurements have been corrected to the Baars flux den- 
sity calibration scale (Baars et al. 1977) using observations of a grid of calibrator sources, 
distributed around the sky, which were observed at approximately two-hour intervals. The 
la error bars include both the estimated measurement and calibration uncertainties. The 
observation and calibration procedures have been described in more detail elsewhere (Aller 
et al. 1985). 

Optical Observations 

In the following, we will discuss three optical data sets. The first was taken at the Boltwood 
Observatory (Stittsville, Ontario, Canada) with a 0.4 m telescope, a SiTe 502A CCD chip, 
and a Johnson- Cousins R filter. The data were collected for ~ 10 days between December 3, 
2002 and January 12, 2003. Relative aperture photometry was performed with an aperture 
of 10 arcsec and "star 1" from Villata et al. (1998) as the comparison star. The background 
was estimated using a concentric annulus with a diameter between 37 and 44 arcsec. We 
did not subtract the light from the host-galaxy. Usually, five two-minute exposures were 
integrated before deriving the photometric value. The typical statistical error on the relative 
photometry is 0.02 mag. The analysis is compromised by two very bright stars near Mrk 421 
that cause a varying level of light to spill into the source and background regions. Based on 
previous optical results on the same source and stars in the field of view, we estimate that 
photometric measurements have a systematic error of 0.08 mag per data point. 

The second set of optical observations were made using a 35 cm Celestron telescope in- 
stalled on the tube of the 60 cm KVA telescope (La Palma, Spain). The observations were 
made with the ST-8 CCD using a standard Kron/Cousins R- filter. The analysis used the 
reduction programs developed by Kari Nilsson (Tuorla Observatory) with the reference stars 
one and two from Villata et al. (1998). 

The final set of optical data were collected using the WIYN 0.9 m telescope at KPNO with 
the S2KB CCD imager using a Harris V-filter. The data were collected from 6 December 
2002 to 15 January 2003 (with, however, a large gap 9 Dec 2002 to 3 Jan 2003). Relative 
aperature photometry was performed using standard IRAF routines with an aperature of 
6 arcsec, sky annulus ranging from 27 to 30 arcsec in diameter, and "star 1" of Villata et 
al. (1998) as the comparison star. Again, we did not subtract light from the host galaxy. 
Typical statistical errors from the photometric fits were smaller than 0.005 mag. Based on 
the spread in magnitude difference between two reference stars, we estimate the uncertainty 
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for each data point to be 0.02 mag for the purpose of determining variability. With regard 
to the absolute flux, due to the presence of host-galaxy light, we expect the values reported 
to contain an undetermined systematic offset of as much as O.f mag. Optical magnitudes 
for all three data sets are converted to absolute fluxes according to Allen (1973). 

X-ray Observations 

We reduced the data from the RXTE Proportional Counter Array with the standard 
RXTE data analysis software. Standard-2 mode PCA data taken with the top layer of 
the operational Proportional Counter Units (PCUs) were analyzed. The number of PCUs 
operational during a pointing varied between 2 and 4. We restricted the spectral analysis 
to the energy range from 4 keV to 15 keV. We excluded data below 4 keV, as the analysis 
of earlier RXTE data showed corrupted behavior (exceptionally high or low count rates of 
individual bins not compatible with the energy resolution of the instrument). We find that 
the data of most pointings are dominated by background above 15 keV. After applying the 
standard screening criteria (including visual inspection of the electron rate) and removing 
abnormal data spikes, the net exposure in each Good Time Interval ranged from 168 sec 
to 9.01 ksec. Light-curves were then extracted with FT00LS v5.3.2. background models 
were generated with the tool pcabackest , based on the RXTE GOF calibration files for a 
"bright" source (more than 40 counts/sec/PCU). Response matrices for the PCA data were 
created with the script pcarsp v. 10.1. The PCU "PCU0" was not excluded for analysis as 
the FTOOLS version gives the proper background model. We assume for all fits a galactic 
neutral hydrogen density of 1.31 xlO 20 cm -2 (Dickey & Lockman 1990) 2 . For each pointing, 
a power law model was fitted over the energy range from 4 keV to 15 keV. 

We complement the data from the pointed RXTE telescopes with data from the RXTE 
All Sky Monitor (ASM) (Levine et al. 1996) taken between December 2, 2002 (UT) and 
January 14, 2003 (UT). We derived fluxes by averaging the "summed band intensities" 
acquired during one day. 



2 http: / / asc.harvard.edu/toolkit / colden.jsp 
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TeV 7-ray Observations 

TeV observations were taken with the Whipple 10 m Cherenkov telescope (Mount Hop- 
kins, AZ) and with the CT1 telescope of the High Energy Gamma- Ray Astronomy (HEGRA) 
collaboration (La Palma, Spain). In the following paragraphs we describe the two data sets. 

The Whipple observations were taken between 4 December 2002 (UT) and 15 January 
2003 (UT). A total of 44 hrs of data were acquired: 32 hrs on the source, and 12 hrs on 
an adjacent field for background estimation purposes. The data were analyzed using the 
standard "Hillas" 2 nd - moment-parameterization technique (Hillas et al. 1985). Standard 
cuts (SuperCuts2000) were used to select 7-ray events and to suppress background cosmic- 
ray events (de la Calle Perez et al. 2003). The fluxes were normalized to the flux from 
the Crab Nebula using a data set of 15 hrs of on-source data and matching background 
observations taken in December 2002 and January 2003 (Punch et al. 1991). Using the 
zenith angle dependence of this Crab data set we account for the zenith angle dependence 
of the 7-ray excess rate by normalizing our measured Mrk421 rate to the Crab rate at 
a corresponding zenith angle. Significances and corresponding error bars were calculated 
using the method of Li & Ma (1983). 

From Monte Carlo simulations, we fold the Crab spectrum with the instrument response 
to obtain the peak energy (energy threshold) for the Whipple 10m data to be consistent 
with the value 400 GeV derived elsewhere (Petry et al. 2002). More detailed descriptions 
of Whipple observing modes and analysis procedures can be found elsewhere (Weekes 1996; 
Punch & Fegan 1991; Reynolds et al. 1993). Details about the Whipple telescope including 
the GRANITE-III camera have been given in Finley et al. (2001). 

A second TeV 7-ray data set was acquired with the HEGRA CT1 telescpe (see Rauterberg 
et al. (1995) for a description of the CT1 instrument) between 3 November 2002 and 12 
December 2002. The telescope was equipped with a 127 pixel camera with a 3 degree 
diameter field of view, and with all-aluminum mirrors giving a total of 10 m 2 reflecting 
surface (Mirzoyan et al. 1994). The analysis used 17 hrs of data with zenith angles between 
12° and 58°. The HEGRA CT1 data were normalized to the Crab flux in a similar way as 
the Whipple data, taking into account the zenith angle dependence of the excess rate. We 
estimate a mean energy threshold for the CT1 data set of approximately 700 GeV. 

The normalization to the steady Crab flux is a convenient way to combine data from 
different instruments to avoid systematic errors resulting from errors in the absolute flux 
calibration of each instrument, and to perform a first-order correction for variations in rate 
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with zenith angle. The drawback of the method is that the different energy thresholds of the 
Whipple and HEGRA observations can introduce a normalization error if the source energy 
spectrum deviates from the Crab energy spectrum. Using previous measurements of the Mrk 
421 TeV spectral index as a function of flux level, we estimate that the Whipple/HEGRA 
normalization error is always smaller than 30%. 

3. Determination of TeV 7-ray Energy Spectra with the Forward Folding 

Method 

The spectral analysis of the Whipple TeV 7-ray data used a different set of gamma- 
hadron separation cuts that minimize the systematic error associated with uncertainties in 
the 7-ray selection efficiency of the cuts while still giving a good sensitivity. The "extended 
zenith angle scaled cuts" (Kosack et al. 2004) select primary 7-rays with an efficiency 
that is largely independent of the zenith angle of the observation and the energy of the 
primary photon. The analysis is based on the Grinnell-ISU (GrISU) package 3 that uses 
the KASCADE airshower simulation code (Kertzman & Sembroski 1994), followed by the 
simulation of the Cherenkov light emitted by the air shower and the simulation of the detector 
response. To calibrate the overall gain of the Whipple 10m telescope in the simulations, we 
compared simulated and observed muon events. Muons show up as bright arcs of Cherenkov 
light in the camera and are useful for calibration because the light per unit arc length is 
nearly constant, regardless of the impact parameter and angle of the muon trajectory. The 
overall gain of the telescope can be found by comparing the distribution of the signal per 
arc length in a simulated set of muon events and in an observed set. We took the simulated 
muon events from a sample of simulated proton and helium showers. Muons are identified 
with a dedicated muon identification algorithm that extracts ~200 muon arcs per 28 min 
data run. We adjusted the overall gain factor in the simulations until they reproduced the 
observed signal per arc length distribution. The overall gain factor agrees to within 15% 
with the value computed from measurements of the mirror reflectivity, photo-multiplier tube 
(PMT) quantum efficiency, and electronic gain. 

We used the forward folding technique to fit the energy spectra. Although earlier TeV 7- 
ray analyses used similar methods, this is the first time that we describe the method in detail. 
For each Cherenkov event that passed the gamma-hadron separation cuts, we computed an 
estimator E of its primary energy, based on the image parameters ' size' S (sum of counts 
in an image) and 1 distance' d (distance of the image centroid from the center of the field of 



3 http://www. physics. utah.edu/gammaray/GrISU 
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view) : 

InE = g(x)+h(d) , (1) 

with x — In S, g(x) = A + Bx + Cx 2 and h(d) = D 1 + a d for d < d and h(d) = D 2 + (3 d 
for d > d . The constants A, B, C, D±, D 2 , a, (3, and d are given in Table 1. The 
first term in the energy estimator reflects the fact that total intensity of an image (size) is 
roughly proportional to the energy (E) of the inducing 7-ray . The second term corrects 
this relationship depending on the distance of the telescope from the shower axis (d is 
proportional to the latter). Using extended zenith angle cuts, the energy estimator gives 
an energy resolution of a (In E) = 0.25. 

We limited the spectral analysis to ON-source data with associated background data sets 
(the so-called ON-OFF data) taken at zenith angles less than 30°. After histogramming the 
energy estimator for both the ON-source and background regions, the background histogram 
was subtracted from the ON-source histogram. Subsequently, an energy spectrum was fitted 
to this "excess histogram" using the forward folding approach (see e.g. Fenimore et al. 
(1982)), making use of a simulated set of 7-ray showers. The simulated set of showers 
consisted of 67,500 showers simulated over an area A M c = 7r(400m) 2 in the energy range 
from 50 GeV to 25.6 TeV over nine energy intervals. The first energy interval went from 
50 to 100 GeV, the second from 100 to 200 GeV and so forth. In the i th energy interval, 
showers were simulated according to a power law distribution: 



dN, 



(0 



MC 



dE 



Ni x (£/lTeV)- rMC (2) 



with T MC = 2.5. Simple integration of Eq. 2 gives the normalization constant AT, as function 
of the lower bound i? m in,i and upper bound E maXti of the i th energy interval and the number 
rii of showers simulated in that energy interval: 

-m x (-Tmc + 1) 
1 F -r MC +i _ F -r MC +i ^ > 

min,i max,i 

We fit the data with two models, a power law model and a power law model with an 
exponential high-energy cutoff: 

|| = iV x(£/lTeV)- r (4) 

and 

^ = N x (E/l TeV)" r x exp (-E/E ) , (5) 

where iVo is the flux normalization at 1 TeV, T is the photon index, and E is the high- 
energy cutoff. For each trial parameter set V (with V = {N ,T} or V = {N ,T,E }), 



another histogram is filled with weighted 7-ray showers that pass the 7-ray selection cuts. 
The weights are computed according to: 



Wi(E;V) 



(6) 
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x (A MC x At) 
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where At is the observation time, and A MC is the area over which showers were simulated. 
The numerator in Eq. 6 gives the model flux for the parameter combination V at energy E. 
The denominator gives the simulated flux. While the weights depend on the true energy of the 
simulated 7-rays , the showers are filled into the histogram according to their reconstructed 
energy. The weighting saves computational time in the fitting procedure, as only the weights 
have to be re-computed for each set of model parameters. 

We performed a search in parameter space until the parameter combination P m i n is found 
that minimizes the x 2 -difference between the observed and simulated histograms. We deter- 
mined the 1 a error region from the condition (Press et al. 1992) 



The best fit model parameters and the associated errors are the main results of a spectral 
analysis. Plotting individual data points in an energy spectrum is well known to be an 
ill-defined problem. A very good discussion can be found in Loredo & Epstein (1989). 
Owing to the finite energy resolution of the telescopes, some information about the true 
energy spectrum is irrevocably lost. We have experimented with "deconvolution methods", 
as for example, the Backus-Gilbert method (Backus & Gilbert 1970; Loredo & Epstein 
1989). Owing to a combination of almost Gaussian distributions of In (E true — E ICC ) and 
E tTue — E Tec and the modest signal-to-noise ratios of the TeV 7-ray energy spectra, we find 
that deconvolution methods improve only very little the effective energy resolution. 

We thus use the simplest method to plot flux estimates, based on the counts in the 
excess histogram (see e.g. Fenimore et al. (1982)). For an energy bin stretching from E\ to 
E 2 , we plot the flux value at the energy E = io^g(E 1 )+iog{E 2 ))/2_ The flux yalue ig giyen by 

scaling the best-fit model according to the observed number of excess counts: 



X\V) < X 2 (Pmin) + 1 • 



(7) 



< ki > 



(8) 



Here, ki is the number of excess counts in the i bin of the signal histogram, and < ki > is 
the sum of weights of the simulated events in the i th bin. 
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Figure 1 shows the Crab spectrum from small zenith angle data (< 30°) acquired be- 
tween 14 September 2002 (UT) and 24 March 2003 (UT), and the time averaged TeV 7-ray 
spectrum of Mrk 421 for the data set of the multiwavelength campaign. Both data sets are 
from the Whipple 10 m telescope. The power law fit to the Crab data gives a flux normaliza- 
tion of 2.42 ± 0.11 x 10" 11 photons cm' 2 s" 1 TeV -1 , and a photon index of T = 2.5 ± 0.1. 
The x 2 is 5.27 for four degrees of freedom. These results should be compared to previous 
results. The Whipple collaboration obtained a flux normalization N = (3.25±0.14±0.60) 
xlO" 11 photons cm" 2 s _1 , and a spectral index T = 2.49±0.06±0.04 (Hillas et al. 1998). 
The HEGRA collaboration published N = (2.83±0.04±0.60) xlO" 11 photons cm" 2 s~\ and 
T = 2.62±0.02±0.05 (Aharonian et al. 2004). Finally, the CAT collaboration obtained N = 
(2.21 ± 0.05 ± 0.60) xlO" 11 photons cm" 2 s -1 , and T = 2.80±0.03 ± 0.06 (Piron et al. 
2001). The first error values are the statistical error values, and the second are the systematic 
error values (both errors on la confidence level). The systematic errors on the absolute flux 
are about 25% and derive from the uncertainty of the energy threshold of the instruments 
combined with the steep spectrum of the Crab Nebula. Our result derived here lies below 
the previously published values of Whipple and HEGRA and above the value published by 
CAT. The discrepancies are comparable to the one-sigma confidence levels. It should be 
noted that the Whipple results published in 1998 were taken with a significantly different 
hardware configuration. Furthermore, we relied here on a energy threshold calibration with 
muons, while Hillas et al. (1998) used the comparison between the detection rate of sim- 
ulated and observed Cosmic Rays to calibrate the energy threshold of the telescope. Each 
calibration method has its own systematic uncertainty, and it is difficult to decide which 
one is more reliable. The discussion shows that one should consider the full systematic error 
when comparing the two Whipple results with each other. The three experiments quote a 
systematic error of ~0.05 on the photon index. Here, our result agrees well with the previous 
Whipple measurements, the HEGRA spectrum is somewhat steeper and the CAT spectrum 
is significantly steeper. The comparison of all the four results shows that systematic errors 
are larger than estimated. In the case of the blazar observations discussed below, the uncer- 
tainty of the absolute energy threshold are not that important, as we are mostly interested 
in relative flux variations. We correct our fluxes with a "throughput factor" derived from 
the Cosmic Ray detection rate measured during each data run to correct for variations in 
the atmospheric conditions. As Cosmic Ray showers are not identical to air showers, we 
estimate that the systematic uncertainty on diurnal fluxes is 10%, and the systematic error 
on diurnal photon indices is 0.1. 



The power law fit to the Mrk 421 data gives a flux normalization of 1.7 ± 0.1 x 10 11 
photons cm" 2 s" 1 TeV" 1 , and a photon index of T = 2.8±0.1. The x 2 is 5.06 for 4 degrees 
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of freedom. The flux and photon index lie in the range of previous observations (Zweerink 
et al. 1997; Aharonian et al. 1999; Krennrich et al. 1999, 2002; Krennrich & Dwek 2003). 
We find that the photon statistics do not allow us to derive meaningful constraints on the 
high-energy cutoff E . 

Using the best-fit parameters from the Crab spectral analysis, we can weight the Monte 
Carlo events by the determined spectrum and compare several simulated distributions of the 
image parameters with the background subtracted distributions for the Crab data (Fig. 2). 
The good agreement between simulated and experimental data verifies that the simulations 
describe the air showers and the detector response adequately. 

4. Results from the Multiwavelength Campaign 

Overview 

Figures 3 and 4 combine all the light-curves measured in December 2002, and January 
2003, respectively. From top to bottom, the figures show the TeV 7-ray data, TeV photon 
indices T, (where dN/dE oc E' r ), RXTE PCU 10 keV flux amplitudes (from the 4-15 keV 
spectral fits), the 4-15 keV photon indices, the RXTE ASM 2-12 keV fluxes, the optical data, 
and the radio data. 

The TeV 7-ray fluxes varied between to ~ 2 times that of the Crab Nebula, with 
slightly higher fluxes observed during the second half of the campaign. We determined 
TeV photon indices on a night to night basis whenever the flux was sufficiently high to 
warrant a spectral analysis. For epochs of low fluxes (December 6 and 7, 2002 (UT) (MJD 
52614-52615), December 8, 9 and 10, 2002 (UT) (MJD 52616-52618), December 14, 15, and 
16, 2002 (UT) (MJD 52622-52624), and January 7 and 8, 2003 (UT) (MJD 52646-52647)) 
we combined the data of several nights to determine an energy spectrum. A x 2 test of 
statistical variability was performed by fitting the entire TeV photon index dataset to a 
constant function. The best fit to the data, with 2a errors, is T = —2.864 ± 0.097, with a x 2 
value of 46.5 for 20 degrees of freedom. The probability of obtaining this value by chance is 
2.2 x 10~ 4 . Seven data points lie outside the 2a confidence region. 

The 4-15 keV RXTE data is shown in Figs. 3 and 4. The X-ray fluxes varied between 
0.2 x 10~ 3 cts cm~ 2 s _1 keV -1 and 4.5 x 10 -3 cts cur 2 s _1 keV -1 . Strong flares were observed 
on December 3, 2002 (MJD 52611), December 5, 2002 (MJD 52613), January 10, 2003 (MJD 
52649), and on January 14, 2003 (MJD 52653). The 4-15 keV X-ray photon indices show a 
large range of values from T = 1.97 to T = 2.90. 
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The values of T >~ 2 indicate that the presumed synchrotron SED peaked at and below 
the 4-15 keV energy range covered by observations. In two cases, the X-ray index varies very 
rapidly: on MJD 52619 it changes by Ar = 0.22 within 1.4 hrs, and on MJD 52650 it 
changes by Ar = 0.14 within 1.6 hrs. 

Over the duration of the campaign, the 4-15 keV photon index changed by Ar « I. 
Synchrotron cooling of a power law distribution of electrons produces a spectral break equal 
to, or smaller than Ar = 0.5 if they high-energy cutoff of the electron spectrum is outside 
the range sampled by the observations. If the power-law index of the accelerated particles 
does not change with time, the detection of spectral variations with Ar > 0.5 thus implies 
that the 4-15 keV fluxes are influenced by the high-energy cutoff of the relativistic electron 
population. Alternatively, flares might be associated with changes of several parameters 
constraining the emitting plasma, as for example the jet magnetic field or the jet beaming 
angle. 

The RXTE rCU and ASM fluxes seem to trace each other, although the sparse sampling 
of the rCU data and the large statistical error bars of the ASM data do not allow us to 
draw definitive conclusions. The optical and radio data show substantial variability by about 
±50% of the mean flux, but no detailed correlation is present. This is not a suprise since 
for the optical data one expects sizeable contributions from the galaxy light and perhaps 
thermal emission from the accretion disk in addition to optical synchrotron. Eurthermore, 
one expects longer variability timescales for the lower energy synchrotron radiation, since 
the lower energy electrons which produce this radiation cool more slowly. 

The X-Ray/TeV 7-ray Flux Correlation 

We studied the X-ray and TeV 7-ray flux correlation with the help of the discrete 
correlation function (DCF) of Edelson et al. (1988). The DCF gives the linear correlation 
coefficient for the two light curves as a function of a time lag between them. The DCF is 
the standard tool used in the case of sparsely sampled data and gives fewer spurious results 
than a traditional correlation function analysis that interpolates between light-curve data 
points. We determined the statistical significance of the correlation coefficients with the help 
of a set of simulated 7-ray light curves, computing for each simulated light curve the DCF 
with the observed X-ray data (Buckley et al. 2005). These light-curves are generated by 
a superposition of triangular shots, with all shots having the same amplitude and the same 
rise and fall time. 
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Figure 5 gives the DCF for the X-ray and TeV 7-ray data sets. For a time lag of 
zero days, we find a DCF value of 0.58 ± 0.12. The simulated data sets show that the 
correlation is significant. For uncorrelated lightcurves consisting of triangular shots with the 
same structure function, we calculate that the chance probability to get a larger DCF value 
at a time lag of zero days is 3.12%. Figure 6 shows the X-ray/TeV 7-ray flux correlation 
for all overlapping observations. The measurements entering this figure are shown by the 
circled data points in Figs. 3 and 4. Each of these overlapping observations are simultaneous 
to within 5 minutes. The correlation shows substantial scatter, with similar X-ray (TeV) 
fluxes sometimes corresponding to TeV (X-ray) fluxes that differ by more than a factor of 
2 from each other. The scatter may be an inherent property of the emission mechanism. 
Alternatively, it may be caused by a short time lag between the flares in the two bands, not 
resolved by the sparse sampling during the campaign. 



Evolution of the X-ray Photon Indices 

Figure 7 shows the correlation between the 10 keV fluxes and the 4-15 keV photon 
indices. A correlation can clearly be recognized in the sense that higher fluxes are accom- 
panied by harder energy spectra. A similar correlation (harder energy spectra for higher 
fluxes) have been reported for other BL Lac-type objects, e.g. Mrk 501 (Kataoka et al. 
1999; Krawczynski et al. 2000), 1ES 1959+650 (Krawczynski et al. 2002), H 1426+428 
(Falcone et al. 2004), and PKS 2155-304 (Kataoka et al. 2000). 

We further scrutinized the temporal evolution of the photon indices with so called X- 
ray "hysteresis" curves (Takahashi et al. 1996; Kirk & Mastichiadis 1999), plotting the 
X-ray photon index as a function of X-ray flux, and indicating the temporal sequence of the 
measurement throughout the evolution of individual flares. In the simplest model whereby 
flares are formed by short lived shocks (e.g. internal shocks from colliding blobs), one 
expects the temporal evolution to be dictated by the interplay of the acceleration, cooling, 
and confinement times. For the Fermi mechanism, the particle energies reached depend 
on the allowed acceleration time. Thus, during the beginning of flares, the X-ray and 7- 
ray energy spectra are expected to harden. Once synchrotron cooling starts to dominate 
the particle energy spectra, the emitted photon energy spectra are expected to soften. We 
studied two X-ray flares, one occurring at the beginning of the campaign (MJD 52612 to 
MJD 52615), and the other at the end (MJD 52651 to MJD 52653). The first flare is shown 
in Fig. 8. This flare coincided with an increase of the TeV 7-ray flux by a factor of 2.4 from 
MJD 52612 to MJD 52613. The X-ray spectrum hardens during the rising phase of the flare 
and softens during the decaying phase. Furthermore, the X-ray spectrum is softer during the 
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falling phase than during the rising phase. The "clockwise" evolution in the T x — F x plane is 
consistent with the expectations from stochastic Fermi acceleration and synchrotron cooling 
as described above. The data of the second flare are shown in Fig. 9. The TeV 7-ray flux 
increased from MJD 52650 to MJD 52651 by a factor of 3.8 and remained roughly constant 
during the following two nights of observations. While the X-ray flux increased from MJD 
52652 to 52653 by a factor of ~ 2.5, the TeV flux measured at the same time as the X-ray 
fluxes did not increase substantially. The general trend is that the spectrum hardens as the 
flux increases, although at MJD 52652.25, the spectrum softens temporarily while the flux 
is still increasing. Unfortunately, our observations did not cover the decaying phase of the 
flare. 



5. Discussion 

The multiwavelength campaign showed Mrk 421 in a level of intermediate activity. During 
the observational campaign, Mrk 421 showed significant flux variability in the radio, optical, 
X-ray and 7-ray bands and significant spectral variability at X-rays and TeV 7-rays . While 
we measured an average TeV 7-ray photon index of T = 2.8, the observations revealed 
evidence for spectral variability on a time scale of days. In particular, the data suggest very 
soft energy spectra with r ~ 4 during the first half of the observation campaign. One of 
the most interesting results from this campaign is that the X-ray and TeV 7-ray fluxes are 
correlated on the ~97% confidence level, but that we find widely different TeV 7-ray fluxes 
for a single X-ray flux and vice versa. The most extreme example is the "orphan X-ray flare" , 
seen on January 13, 2003 (MJD 52653). The loose X-ray/TeV correlation may suggest that 
the model parameters (e.g. the volume of the emission zone) change with time, or that the 
commonly made assumption of a single synchrotron self-Compton (SSC) emission zone is 
an over-simplification. Our previous observations of Mrk 421 had already shown a rather 
loose X-Ray/TeV 7-ray correlation (Blazejowski et al. 2005) and the same applied for 1ES 
1959+650 (Krawczynski et al. 2004). In the case of Mrk 501, a rather tight correlation has 
been reported (Krawczynski et al. 2000, 2002). 

The analysis of the correlation between the X-ray flux and photon index during a flare 
indicated a "clockwise" hysteresis. For Mrk 421, Takahashi et al. (1996) also reported 
clockwise evolution. However, Takahashi et al. (2000) reported evidence for both, clockwise 
evolution during some flares and anti-clockwise evolution during other flares. If the SSC 
model indeed applies, these results may imply that the characteristic times scales of the most 
important processes (acceleration time, radiative cooling time, escape time) change from flare 
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to flare and thus yield the different observed hysteresis behaviours. Recently, Sokolov et al. 
(2004) emphasized that the geometry of the emitting region and its orientation relative to 
the line of sight influences the observed flux and spectral evolution and might thus further 
complicate the interpretation of the results. 

The X-ray and TeV 7-ray emission from Mrk 421 data have been modeled with syn- 
chrotron Self-Compton codes by many groups, see e.g. (Inoue & Takahara 1996; Bednarek 
1997; Bednarek & Protheroe 1999; Bottcher, Mause & Schlickeiser 1997; Mastichiadis & 
Kirk 1997; Tanahita et al. 2001; Krawczynski et al. 2001; Konopelko et al. 2003; Kino et al. 
2002; Blazejowski et al. 2005). A crucial model parameter is the jet Doppler factor 8j. The 
published models with Doppler factors 5j of 20 or lower generally predict TeV energy spectra 
that are softer than the observed ones, especially if a correction for extragalactic absorp- 
tion would be applied which steepens the energy spectra considerably. Models with Doppler 
factors <5j on the order of 50 give satisfactory model fits to both the X-ray (synchrotron) 
and the TeV (Inverse Compton) emission (see the self-consistent modeling of Krawczynski 
et al. (2001); Konopelko et al. (2003) and the discussions by Tavecchio (2004); Piner 
& Edwards (2005)). Piner & Edwards (2005) observed the Mrk 421 parsec-scale radio 
jet with the Very Large Baseline Array (VLBA). Remarkably, they find apparent pattern 
speeds of only ~ 0.1c. As discussed by the authors, the highly relativistic motion inferred 
from TeV observations can be reconciled with the modestly relativistic flow calculated from 
VLBA observations by postulating that the jet slows down between the sub-parsec (TeV) 
and parsec (VLBA) regimes. It may be possible to describe the multi-wavelength data with a 
synchrotron-Compton model and lower Doppler factors by invoking additional seed photons. 
Two new model variants that combine ingredients of SSC and external Compton models have 
been proposed by Georganopoulos & Kazanas (2003) and by Ghisellini et al. (2005). While 
the first authors assume that downstream emission regions provide seed photons, the second 
authors speculate that the jet is composed of a fast spine with a slow-moving envelope. In 
this model, the fast spine emits the X-ray and 7-ray radiation. Modeling of the data with 
these two inhomogeneous models is outside the scope of this paper. 

In this discussion we do not want to embark on comprehensive modeling of the data 
from the entire campaign. Our main aim is to draw the attention of the reader to a single 
remarkable fact: while it is difficult to model the data with Doppler factors on the order of 
20 and lower, much higher Doppler factors cannot be excluded right away. Fig. 10 shows 
two synchrotron self-Compton models based on the simple snapshot code of Krawczynski et 
al. (2004). The code assumes a single spherical emission volume of radius R relativistically 
approaching the observer with a jet Doppler factor <5j. The emission volume is homogeneously 
filled with a tangled magnetic field of strength B and a non-thermal electron population. 
The electron energy spectrum follows dN/d'y oc 7~ p with p = 2 for electron Lorentz factor 
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7 between 7 min and 7b and p = 3 for Lorentz factors between 7b and 7 max - The code models 
extragalactic absorption owing to the 7x e v + 7cib —> e + e~ pair-production processes of TeV 
photons on photons from the Cosmic Infrared Background (CIB) using the CIB model of 
Kneiske, Mannheim & Hartmann (2002). 

We discuss two models. We show the first model for illustrative purposes only. It uses the 
"conventional" model parameters (5j = 50) derived from the time-dependent self-consistent 
analysis of a different but similar data set (Krawczynski et al. 2001). The second model uses 
a very high Doppler factor (5j = 1000). All the model parameters are given in Table 2. For 
both models, we assured that the model parameters were chosen self-consistently. Causality 
arguments require that the radius R satisfies R < 5jcAT ^ s = 2.7 x 10 15 cm for 5j = 50 and 
R < 5.4 x 10 16 cm for 5 i = 1000 for a flux variability time scale of AT obs = 30 min. Note 
that the flux variability time scale sets a lower limit on R but no upper limit, if the flux 
variability time scale is not dominated by light travel time effects but by other effects (e.g. 
by the stability of a strong shock front). We checked that the SED (i) fits the X-ray and TeV 
7-ray data, and (ii) is consistent with the expected spectral shape owing to radiative cooling. 
In the first model, the latter self-consistency is assured by our previous self-consistent time- 
dependent modeling. In the second model, we construct an electron spectrum based on the 
general results for electrons suffering synchrotron and Inverse Compton losses (Syrovatskii 
1959; Kardashev 1962; Ginzburg & Syrovatskii 1964; Pacholczyk 1970; Inoue & Takahara 
1996)). We assume that the electron energy spectrum breaks at 7b = 1.8 x 10 3 . The 
laboratory-frame synchrotron cooling time for electrons at the break is t s = [| cr-p c 5j 
g7r ra e 2 7b ] _1 ~ 28mm (a T is the Thomson cross section and m e is the electron mass). An 
electron spectrum as the one used here could result from the radiative cooling of a p =2 
electron energy spectrum that extends from 7 min to 7 max! resulting in a spectrum with p = 2 
and p = 3 below and above 7b, respectively. 

As shown in Fig. 10, both models fail to predict the observed radio fluxes as a conse- 
quence of synchrotron self-absorption. We would like to emphasize that we do not regard 
the discrepancy as a shortcoming of the model. Electrons producing the radio emission cool 
on much longer time scales, and the radio emission will be dominated by an accumulation 
of downstream plasma which stopped contributing to the X-ray and TeV emission long time 
ago. We could model the radio emission with another emission component (see e.g. Blaze- 
jowski et al. (2005)). However, doing so is arbitrary: for a small number of data points we 
would add an additional model component with many free model parameters. 

In Table 2 we list for both models the magnetic field energy density ub, the energy 
density in electrons u e , the ratio r = u e /uB, and the kinetic luminosity = tiR 2 cT 2 (u c + ub) 
for r = <5j (Begelman et al. 1994). The model with a low Doppler factor is closer to 
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equipartition between magnetic field and particles. The kinetic luminosity is similar for the 
two models, with a high radiative efficiency of the low-<5j model making up for the stronger 
boosting of the high-<5j model. 

If taken seriously, the model with 5j = 1000 would imply that the X-ray and TeV 7- 
ray emission is produced by an ultra-relativistic particle dominated wind, very close to the 
supermassive black hole. The fact that seven blazars have been detected at TeV energies 
seems to argue against extremely relativistic outflows with bulk Lorentz factors T on the 
order of 1000, as isotropic emission would be beamed into an opening angle of T _1 and 
would make the observation of the emission unlikely. However, the argument only applies 
if the jet opening angle is equal or smaller than r -1 . Having a larger jet opening angle 
would require a higher total jet-luminosity as some jet segments would not contribute to 
the observed emission. However, the jet-luminosities listed in the table are several orders 
of magnitude below the Eddington luminosity of a ~ 10 8 ' 4 solar mass black hole that is 
suspected to be at the core of Mrk 421 (Barth, Ho & Sargent 2003; Falomo, Kotilainen & 
Treves 2002). 
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Fig. 1. — Whipple TeV spectrum of the Crab Nebula and Mrk 421. The dashed and dotted 
lines give the results of power law fits for the Crab and Mrk 421, respectively. 
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Fig. 2.— Hillas parameter distributions for Whipple 2002-2003 Crab ON/OFF data and 
simulated data. Histograms show simulated data, while data points show Crab data. 
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Fig. 3. — Multi-wavelength data from December, 2002. The 7-ray data points show per- 
night averages, in Crab units. Starred data points signify data taken from the HEGRA CT1 
telescope, while crosses denote Whipple 10m data. The error bars on the RXTE PCA data 
are not shown as they are smaller than the symbol size, and have units of 10 -10 erg cm -2 
s _1 at 10 keV. The circled X-ray and 7-ray data points overlapped or were taken less than 
5 min apart. The TeV 7-ray and RXTE photon indices show T, where dN/dE oc E~ r . the 
ASM data are given in mCrab. In the optical band, open circles show the WIYN V band 
data, crosses show the Boltwood R band data, and 'x' denotes La Palma R band data. All of 
the optical data are in relative magnitude units. The open (filled) circles in the radio band 
show measurements that overlapped or were taken within 5 min of a TeV 7-ray observation 
(X-ray and TeV 7-ray observation). The radio data are given in Janskys. 
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Fig. 4.- 



Same as Fig. 3 for the data from January 2003. 
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Fig. 5. — Discrete correlation function of the complete X-ray and 7-ray data set. A positive 
time lag means the 7-ray flux precedes the X-ray flux. 
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Fig. 6. — Plot of the TeV 7-ray versus X-ray flux correlation for measurements for all 
overlapping data sets. 
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Fig. 7. — Correlation of the 10 keV X-ray flux and the3-14 keV photon index (both: RXTE 
PC A data). 
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Fig. 8. — X-ray power law photon index versus X-ray flux for the MJD 52612-52615 X-ray 
flare. 
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Fig. 9. — X-ray power law photon index versus X-ray flux for the MJD 52651-52653 X-ray 
flare. 
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Fig. 10. — Mrk 421 Spectral Energy Distributions measured during the campaign. The 
data points show the radio to 7-ray data. The radio and optical spectra show the average 
fluxes observed during the campaign. At X-rays we show three spectra, one from the RXTE 
pointing with the highest flux observed during the campaign, one from the pointing with the 
lowest flux, and one spectrum at intermediate flux levels. The intermediate X-ray spectrum 
and the gamma-ray spectrum were determined using for both only the data taken during 
nights with simultaneous X-ray and gamma-ray observations. The solid curved lines show, for 
comparison, a low-flux and a high-flux energy spectrum measured with BeppoSAX ( "Beppo" 
Satellite per Astronomia X) during the 1998 flaring period (Fossati et al. 2000). the long 
dashed lines show the results from a simple Synchrotron Self-Compton model with 5 = 1000, 
while the short dashed lines show the results with 5 = 50. 
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Parameter Value 



A 


-7.05 


B 


1.30 


C 


-0.034 


D l 


0.057 


a 


-0.20 


D 2 


-1.96 


P 


2.44 


d 


0.75° 



Table 1: Parameters used for the TeV 7-ray energy estimator. 
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8 B (G) R (m) log7 min log7 max log7 b U B U part -^f L k 



50 0.2 1.05 x 10 13 3.3 5.3 4.8 1.92 x 10~ 3 8 x 10~ 2 41.5 2.12 

1000 0.5 1.7 x 10 10 2.8 4.4 3.26 9.95 x 10~ 3 400 4.02 x 10 4 10.9 

Table 2: Parameters for 2 synchrotron self-Compton models. Here, 8 is the relativistic 
Doppler factor, B is the magnetic field, in Gauss, R is the size of the emission region, in 
meters, r y m in,'lmax and 7& are the starting, ending, and break Lorentz factors for the primary 
electron energy spectrum, Ub and U part are energy densities in units of erg cm" 3 , and is 
the minimum kinetic luminosity (defined in section 5) in units of 10 43 erg s _1 . 



